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Abstract 

The coherent vortex extraction method, a wavelet technique for extracting coherent vortices 
out of turbulent flows, is applied to simulations of resistive drift-wave turbulence in magnetized 
plasma (Hasegawa-Wakatani system). The aim is to retain only the essential degrees of freedom, 
responsible for the transport. It is shown that the radial density flux is carried by these coherent 
modes. In the quasi-hydrodynamic regime, coherent vortices exhibit depletion of the polarization- 
drift nonlinearity and vorticity strongly dominates strain, in contrast to the quasi-adiabatic regime. 
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I. INTRODUCTION AND GOVERNING EQUATIONS 



One important issue in fusion research is the understanding and control of turbulent ra- 
dial flux of particles and heat in magnetized plasmas, in order to improve the confinement 
properties of fusion devices^. Indeed turbulence enhances the radial diffusion dramatically 
compared to neo-classical estimations. A long standing question has been^""-: what is the 
role of coherent structures in this radial transport? The answer to this question requires 
extracting and characterizing coherent structures. A particularly appropriate framework to 
identify coherent structures is the wavelet representation, where wavelets are basis functions 
well localized in both physical and Fourier spaced. It has already been used to identify 
coherent structures in fluid turbulence and to distinguish them from background incoherent 
noise^. These methods have recently been applied to experimental signals of ion density 
in the tokamak scrape-off layer~, separating coherent bursts from incoherent noise. In the 
present work these methods are applied to assess the role of coherent vorticity structures 
in anomalous radial transport in two-dimensional numerical simulations of drift-wave tur- 
bulence. Drift waves are now generally considered to play a key role in the dynamics and 



transport properties of tokamak edge turbulence (e.g. [10|] and references therein). At the 
edge, the plasma temperature is low and the collision rate relatively large, therefore the 
resistivity is potentially important. The Hasegawa-Wakatani mode\^>^ is a two-field model 
which includes the main features of turbulent transport by resistive drift waves. 

In the present work the two-dimensional slab geometry-version of this model is chosen as 
a paradigm for drift-wave turbulence in the plasma-edge region. In dimensionless form the 
Hasegawa-Wakatani model reads^^ 

^ - ^V^) n + K^ + cin-<P) = [n, 0] , (1) 

^ - z/V^^ + c(n - 0) = [V20, 0] , (2) 

with n the plasma density fluctuation and the electrostatic potential fluctuation. D 

and u are the cross-field diffusion of plasma density fluctuations and kinematic viscosity, 

respectively. The Poisson brackets are defined as 

da db da db . , 

dx dy dy dx' 

We identify the x-coordinate with the radial direction and the y-coordinate with the poloidal 
direction. The equilibrium density no is non- uniform, with a density gradient dno/dx in the 



negative x-direction, such that the equihbrium density scale L„ = no / {duQ / dx) is constant 
and the value of k is one. The plasma density fluctuations n are normalized by no, therefore 
n/riQ — )■ n, the electrostatic potential is normalized as ecp/T^ -> 0, the space as x/ps x 
and the time as Udt — )■ t, where e is the electron charge, Te the electron temperature, Ud 
the ion cyclotron frequency and ps = {miTeY^'^/ (eB) is the ion integral lengthscale. B is the 
strength of the equilibrium magnetic field in the 2;-direction and is the ion mass. The key 
parameter in this model is the adiabaticity c, which represents the strength of the parallel 
electron resistivity. It is defined as 



(4) 



with the effective parallel wavenumber and 1] the electron resistivity. 
The vorticity u is related to the electrostatic potential (p by 

= (5) 

Note that for c = 0, equation ([1]) corresponds to the advection-diffusion of a passive scalar in 
the presence of a (unity) mean scalar gradient in the x-direction. Equation ([2]) corresponds in 
this case to the vorticity equation. For c — )■ oo the Hasegawa-MimaM one field approximation 
is approached,— which ignores all resistive effects. For c — )■ we recover the hydrodynamic 
limit, which is less relevant to describe edge fusion-plasma. Here two cases will be considered: 
a quasi-adiabatic case with c = 0.7, and a quasi-hydro dynamic case with c = 0.01. The case 
c = 0.7 is generally considered to be the most relevant for tokamak- research and has been 
investigated in several other works (e.g. {4, 15]). Both cases differ from the fluid-dynamical 



case in that the velocity field is forced through the interaction term c{n — (p). The influence 
of this term on the density field can however be considered to be negligible in the quasi- 
hydrodynamic case^. 

The quantity of interest, the radial particle density flux, is the correlation between the 
radial velocity Ur = —dcp/dy and the particle density, 

Fj. = {nur) , (6) 

where the brackets denote an average over both time and space. The question we address in 
this paper is how coherent structures contribute to this flux. To investigate this, direct nu- 
merical simulations of the Hasegawa-Wakatani system are performed on a periodic domain 
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FIG. 1: One realization of the vorticity field for the quasi- hydrodynamic case (left) and for the 
quasi-adiabatic case (right). The abscissa corresponds to the radial position. The ordinate indicates 
the poloidal position. Both range from to 64 ps- The white frames indicate the dipoles we have 
selected in both cases. 

discretized with = 512^ gridpoints. The length of the domain is 64 p^. A finite difference 
method is used in which the nonlinear terms are computed using a method developed by 
Arakawar^. The time stepping is performed using a predictor-corrector scheme. The plasma 
density diffusion D and viscosity v are set to 0.01 in normalized units. Computations are 
performed up to t = 612. At t ~ 100 the kinetic energy saturates and a statistically station- 
ary state is reached, independent of the (random) initial conditions. Typical realizations 
of the vorticity field are shown in figure [T|, where one observes coherent structures for both 
cases. In each case we select a dipolar structure that we indicate by a white frame. The 
quasi-hydrodynamic case exhibits coherent vortices of very different sizes and intensities, in 
contrast to the quasi-adiabatic case where the coherent structures are more similar in size 
and intensity. 



II. COHERENT VORTEX EXTRACTION (CVE) 
A. Method 



Definitions and details on the orthogonal wavelet transform and its extension to higher 
dimensions can be found, e.g., in 0il7|. In the following we fix the notation for the orthogonal 
wavelet decomposition of a two-dimensional scalar valued field. The wavelet transform 
unfolds the field into scales, positions and directions using a set of dilated, translated and 
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rotated functions, called wavelets. Each wavelet is well-localized in space, oscillating {i.e., 
it has at least a vanishing mean, or better its first m moments vanish), and smooth {i.e., its 
Fourier transform exhibits fast decay for wavenumbers tending to infinity). We here apply 
the coherent vortex extraction (CVE) algorithnt^ii^ using orthogonal wavelets. In dimension 
two, orthogonal wavelets span three directions (horizontal, vertical and diagonal), due to 
the tensor product construction. To go from one scale to the next, wavelets are dilated 
by a factor two and the translation step doubles accordingly. Wavelet coefficients are thus 
represented on a dyadic grici^. 

We apply the CVE algorithm to the vorticity fields u of both the quasi-hydrodynamic and 
the quasi-adiabatic regime. The extraction is performed from the vorticity since enstrophy 
is an inviscid invariant in the hydrodynamic limit. Moreover, vorticity is Galilean invariant 
in contrast to velocity and streamfunction. We consider the quasi-stationary state of the 
simulations, i.e., when a saturated regime is reached, and we decompose the vorticity field, 
given at resolution = 2^"^, into an orthogonal wavelet series 



where the multi-index A = {j,ix,iy,d) denotes the scale j the position i = 
{ix,iy) and the three directions d = 1,2,3, corresponding to horizontal, verti- 
cal and diagonal wavelets respectively. The corresponding index set A is A = 
{A = {j, ix, iy, d),j = 0, J — I; ix, iy = 0...2^ — 1, d = 1,2,3} . Due to orthogonality the 
wavelet coefficients are given by u\ = {u,ipx), where (■,■) denotes the L^-inner product 
defined as {f,g) = j f{x,y)g{x,y)dxdy. The wavelet coefficients measure fiuctuations of 
u at scale around the position i, in one of the three directions d. Here a Coifman 30 
wavelet is used, which is orthogonal and has 10 vanishing moments^^. {J x^ip{x)dx = for 



The CVE algorithm can be summarized in the following three step procedure: 

• Decomposition: compute the wavelet coefficients u\ using the fast wavelet transform^. 

• Thresholding: apply the thresholding function to the wavelet coefficients uj\, thus 
discarding the coefficients with absolute values smaller than the threshold e. 

• Reconstruction: reconstruct the coherent vorticity field uc from the thresholded 




(7) 



n = 0,...9). 
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TABLE I: Compression rate (% of coefficients retained), retained energy E = ^ (^ip ^'^'P): 
enstrophy Z = ^ ('^^)' ^'^'^ radial flux T^, after applying the CVE filter to the vorticity field of 
the quasi-hydrodynamic and quasi-adiabatic 2D drift-wave turbulence simulations. 





Compr. ( 


;%) E (%) z (%) r, (%) 


Quasi-hydrodynamic (c=0.01) 


1.3 


99.9 97 99 


Quasi-adiabatic (c=0.7) 


1.8 


99.0 93 98 



wavelet coefficients using the fast inverse wavelet transform. 

The incoherent vorticity field is obtained by simple subtraction, i.e., cj/ = a; — uc- 
The thresholding function is given by 

Pe{a) = , (8) 

where e denotes the threshold, 

e = VAZlnN, (9) 

where Z = u) is the enstrophy (which corresponds to half of the variance of the vorticity 
fiuctuations) and the resolution. This threshold value allows for optimal denoising in a 
minmax sense, assuming the noise to be additive, Gaussian and white^. 

In summary, this decomposition yields u = uq + w/. Due to orthogonality we have 
{uc^oji) = and hence it follows that enstrophy is conserved, i.e., Z = Zc + Zj. Let us 
mention that the computational cost of the Fast Wavelet Transform (FWT) is of 0{N)^. 

B. Compression rates 

The results of the extraction are displayed in table [H The compression rate is in both 
cases very significant: for the quasi-hydrodynamic case, 1.3% of the modes retain more than 
99.9% of the energy and 97% of the enstrophy. For the quasi-adiabatic case, 1.8% of the 
modes retain 99.0% of the energy and 93% of the enstrophy. The contribution of the coherent 
vorticity to the radial fiux is also given in table HI The coherent modes, which contain most 
of the energy and enstrophy, are responsible for 99% of the radial particle density fiux in 
the quasi-hydrodynamic case, and for 98% of F^ in the quasi-adiabatic case. In other words, 
F^ is almost exclusively carried by the coherent structures. 
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FIG. 2: Top: PDF of the vorticity. Bottom: Fourier spectrum of the enstrophy versus wavenumber. 
Left: quasi-hydrodynamic case. Right: quasi-adiabatic case. Dashed hne: total field, solid line: 
coherent part, dotted line: incoherent part. Note that the coherent contribution (solid) superposes 
the total field (dashed), which is thus hidden under the solid line in all four figures. The straight 
lines indicating power laws are plotted for reference. 

C. Wavenumber spectra and probability density functions 

Spectra and probability density functions (PDF), averaged over 512 realizations during 
the time interval 100 < t < 612, are shown in Figure[2]for the total, coherent and incoherent 
vorticity. The PDF of the total and coherent quasi-hydrodynamic vorticity is far from Gaus- 
sian and slightly skewed, while the quasi-adiabatic vorticity is much closer to Gaussianity. 
In both cases, the variance of the incoherent part is much smaller than the variance of the 
coherent part, which has the same PDF as the total. For the quasi-hydrodynamic case, the 
coherent part retains 97% of the variance of the vorticity fluctuations and therefore also 
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FIG. 3: (Color online) Scatter-plot of vorticity against electrostatic potential for the coherent part 
(top) and incoherent part (bottom). Left quasi- hydrodynamic case, right quasi-adiabatic case. The 
light grey (red online) dots correspond to the total field, the dark grey (blue online) dots to the 
dipoles we have selected in Fig. 1. 

97% of the total enstrophy Z, with Z = 1.4. For the quasi-adiabatic case, the coherent part 
retains 93% of the variance of the vorticity fluctuations and hence 93% Z, with Z = 3.4. A 
similar result is observed in the enstrophy spectrum computed from the Fourier transform of 
the vorticity field, averaged over wavenumber shells of radius |A;|, the wavenumber. The total 
and coherent enstrophy are the same all over the inertial range and at the highest wavenum- 
bers, in the dissipation range, the incoherent part contributes to the spectral enstrophy 
density. Both coherent and incoherent contributions are spread all over the spectral range, 
but they present different spectral slopes in the inertial range and therefore different spatial 
correlations. From the integral wavenumber to the dissipation wavenumber, a negative slope 
for the coherent contribution, corresponding to long range spatial correlations, is observed. 
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The incoherent part shows a positive slope with a power-law dependence close to in the 
inertial range. This corresponds to an equipartition of kinetic energy in two dimensions. A 
similar result was obtained in three-dimensional isotropic Navier-Stokes turbulence^. 

D. Scatter-plots 

We show in figure [3] scatter-plots of the vorticity versus the electrostatic potential corre- 
sponding to the fields in figure [H Both the total part and the incoherent part are shown. 
Since the coherent part is almost identical to the total part, it has been omitted. Also shown, 
superposed on the same figures, is the scatter-plot corresponding to the zoom on the dipo- 
lar structures indicated by a white frame in figure [H In the freely decaying hydrodynamic 
case, c = 0, Joyce and MontgomeryiS showed that a functional relation (j){uj) = asinh.{Pu) 
should be expected, corresponding to a final state of decay depleted from nonlinearity. The 
parameters a and (3 are Lagrangian multipliers, necessary for maximizing the entropy under 
constraints. The value 1/(3 can be associated with a (negative) temperature^-. Depletion 
from nonlinearity corresponds to steady solutions of the Euler equation, [u, 0] = 0, implied 
by the existence of a functional relation <f){u). Indeed drift- wave turbulence contains an 
internal instability which prevents the fiow from decaying. This forcing is present in both 
cases considered here and a sinh-Poisson relation cannot be expected a priori for the global 
fiows. Moreover, the two-field model [equations ([T]) and ([2])] contains two nonlinearities, first 
the polarization-drift nonlinearity in the vorticity equation, second the E x B nonlinearity 
in the density equation. The latter disappears in the adiabatic limit as n and are in phase, 
which corresponds to a linear functional relationship. In figure [3l a local depletion of the 
polarization-drift nonlinearity is seen for the quasi-hydro dynamic case. The scatter-plot of 
(j) — oj, corresponding to the dipolar structure, that is indicated by a white frame in figure 
[U (left), is close to a sinh-Poisson relation (solid black curve) in spite of the presence of 
the forcing term. In the quasi-adiabatic case the dipolar structure, that is indicated by a 
white frame in figure [1] (right), does not exhibit such a functional relation. In the incoher- 
ent parts (Figure |3l bottom) no functional relation can be distinguished, which confirms 
that the incoherent part does not contain any structure, for both quasi-hydro dynamic and 
quasi-adiabatic cases. 
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E. Strain versus vorticity 



A question is now how to quantitatively distinguish between the structures in both cases. 
Intuitively it can be inferred that different regions of high vorticity in the quasi-adiabatic 
case involve strong mutual shearing which strongly limits their lifetime and the chance to 
reach a functional relation (f){u)). Koniges et al. - determined the lifetime of individual eddies 
compared to the eddy-turnover time Tover, i-e. the time it takes for a fluid element in an eddy 
to make a 27r rotation. They estimated the lifetime of the quasi-hydrodynamic eddies to be 
approximately 10 Tg^er, and the lifetime of the adiabatic eddies (for c = 2.0) approximately 
Toyer- As mentioned in their paper, this measure is quite subjective and very time-consuming, 
especially if a full PDF of the lifetimes is to be obtained. Here we propose a simpler way to 
distinguish the coherent structures for the different regimes. 

In fluid turbulence the Weiss criterion is a local measure of the strain compared to 
the vorticity for a 2D velocity field. The Weiss field is defined as: 

Q = l{a'- uj^) , (10) 

with 

\dx dy J \dy dx J 
u and V are two orthogonal components of the velocity vector. The Weiss criterion was 
proposed to identify coherent structures, but it may lead to ambiguous results because 
the underlying assumption that the velocity gradient varies slowly with respect to the vor- 
ticity gradient is not generally valicl^i. We here apply the same criterion to drift-wave 
turbulence22i"— but not to identify coherent structures (this being done by the CVE method), 
but to distinguish between the quasi-hydrodynamic and quasi-adiabatic cases. 

The PDF of the Weiss field (FigJl]) reveals that it is its skewness that differentiates best 
the two fields. Indeed, it is more skewed towards negative Q for the quasi-hydrodynamic 
case than for the quasi-adiabatic case: the skewness is —11 for the former, compared to —2 
for the latter. The PDF shows thus that in the quasi-hydrodynamic case the probability to 
find rotationally dominated regions is larger, and the rotation exhibits much larger values, 
than in the quasi-adiabatic case. The variance of Q is comparable for the two cases (5 and 4 
for the quasi-hydrodynamic case and the quasi-adiabatic case, respectively). The skewness 
of the Weiss field Q appears to be a good quantitative measure to distinguish between the 
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FIG. 4: PDF of the Weiss field Q for the quasi-hydrodynamic and quasi-adiabatic velocity fields. 

two cases studied in the present work. In further studies it can be investigated, whether this 
measure can be used to identify coherence in different types of turbulent flows. 

III. CONCLUSION AND PERSPECTIVES 

In conclusion, we have applied the Coherent Vortex Extraction method to dissipative 
drift-wave turbulence. The results show that we can identify the essential degrees of freedom 
(less than 2%) responsible for the nonlinear dynamics and transport. The coherent modes 
contain almost all the energy and enstrophy and contribute to more than 98% of the radial 
flux. 

Evaluating the scatter-plot of the vorticity versus the electrostatic potential, it is shown 
that the coherent structures in the quasi-hydrodynamic case are close to a state of local 
depletion of polarization-drift nonlinearity. In contrast, this is not the case for the quasi- 
adiabatic regime, where nonlinearity remains active and no sm/i-functional relation between 
vorticity and electrostatic potential is observed. This depletion of nonlinearity in the quasi- 
hydrodynamic regime may explain the failure of the quasi-linear estimate of the radial flux-. 
The skewness of the Weiss field yields a quantitative measure for the difference in nonlinear 
behavior of the coherent structures between the quasi-hydrodynamic and quasi-adiabatic 
cases. 

The wavelet transforms, or the proper orthogonal decomposition (POD), may become 
very useful to denoise particle-in-cell simulations of plasma turbulence^^. A comparison of 
the performance of the POD and CVE method is currently undertaken and will be reported 
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in a future paper. 
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